/**
  CNOK project, Anyang Normal University, IMP-CAS
  \brief Calculate inclusive Rs=sigma_exp/sigma_th for -N reactions.
  \author SUN Yazhou, asia.rabbit@163.com
  \since 2022/05/23
  \date Last modified: 2022/05/23 by SUN Yazhou
  \copyright 2020-2023 SUN Yazhou
  \copyright CNOK project, Anyang Normal University, IMP-CAS
*/

#include <iostream>
#include "TARs.h"
#include "TAMOMDIS.h"
#include "TAParaManager.h"
#include "TAException.h"
#include "TANucleus.h"

using std::cout;
using std::endl;

#define O16
static const int FNNPAR = 1; // 0: Horiuchi; 1: LenziRay
static const bool PAULI = false; // 0: off; 1: on
// static const double EK = 2100.; // incident energy: MeV/nucleon

int main(int argc, char const *argv[]){
  // n, sp, sn, cs, dcs
  TANucleus o13 = TANucleus(13, 1.5116, 16.8705, 21.56, 3.34);
  o13.SetBaseDir("config/O/O13");
  o13.SetOrbit({"0p1/2_1+", "0p3/2_1+"});
  o13.SetEx({0., 0.});
  o13.SetC2S({0.5291, 0.0923});

#ifdef O13
  // path to the yaml files
  const char *basedir = "config/O/O13";
  const double A = 13., Sp = 1.5116, Sn = 16.8705;
  const double csExp = 21.56, dcsExp = 3.34;
  const int n = 2;
  const char *orbit[n] = {"0p1/2_1+", "0p3/2_1+"};
  const double Ex[n] = {0., 0.};
  const double c2s[n] = {0.5291, 0.0923};
#endif

#ifdef O14
  const char *basedir = "config/O/O14";
  const double A = 14., Sp = 4.6267, Sn = 23.1790;
  const double csExp = 38.65, dcsExp = 3.01;
  const int n = 1;
  const char *orbit[n] = {"0p1/2_1/2-"};
  const double Ex[n] = {0.};
  const double c2s[n] = {1.6102};
#endif

#ifdef O15
  const char *basedir = "config/O/O15";
  const double A = 15., Sp = 7.2968, Sn = 13.2235;
  const double csExp = 65.38, dcsExp = 6.06;
  const int n = 14;
  const char *orbit[n] = {"0p1/2_1+_1", "0p1/2_1+_2", "0p1/2_1+_3", "0p1/2_0+",
  	"0p3/2_1+_1", "0p3/2_1+_2", "0p3/2_1+_3", "0p3/2_2+", "0d5/2_2-", "0d5/2_3-",
  	"1s1/2_0-", "1s1/2_1-", "0d3/2_2-", "0d3/2_1-"};
  const double Ex[n] = {0., 3.9481, 6.2035, 2.3128, 0., 3.9481, 6.2035, 7.0291,
  	5.1059, 5.8343, 4.9151, 5.6914, 5.1059, 5.6914};
  const double c2s[n] = {0.7972, 0.5401, 0.005, 0.4052, 0.2998, 0.3678, 0.016,
  	1.15, 0.0358, 0.0257, 0.0001, 0.0074, 0., 0.0004};
#endif

#ifdef O16
  const char *basedir = "config/O/O16";
  const double A = 16., Sp = 12.1274, Sn = 15.6639;
  const double csExp = 60.11, dcsExp = 3.38;
  const int n = 10;
  const char *orbit[n] = {"0p1/2_1/2-", "0p3/2_3/2-_1", "0p3/2_3/2-_2",
  	"0d5/2_5/2+_1", "0d5/2_5/2+_2", "0d3/2_3/2+_1", "0d3/2_3/2+_2",
  	 "1s1/2_1/2+_1", "1s1/2_1/2+_2", "1s1/2_1/2+_3"};
  const double Ex[n] = {0., 6.3238, 9.925, 5.2702, 7.1551, 7.3008,
  	8.5714, 5.2988, 8.3126, 9.0497};
  const double c2s[n] = {1.8049, 3.7203, 0.0251, 0.0961, 0.0141, 0.0036,
  	0.0235, 0.0081, 0.0095, 0.0001};
#endif

#ifdef O17
  const char *basedir = "config/O/O17";
  const double A = 17., Sp = 13.7816, Sn = 4.1431;
  const double csExp = 21.4, dcsExp = 3.8;
  const int n = 5;
  const char *orbit[n] = {"0p1/2_2-", "0p1/2_3-", "0p3/2_2-", "0p3/2_1-", "0p3/2_3-"};
  const double Ex[n] = {0., 0.2982, 0., 0.3973, 0.2982};
  const double c2s[n] = {0.7225, 1.0365, 0.0246, 0.007, 0.0173};
#endif

#ifdef O18
  const char *basedir = "config/O/O18";
  const double A = 18., Sp = 15.9419, Sn = 8.0454;
  const double csExp = 29.6, dcsExp = 2.2;
  const int n = 9;
  const char *orbit[n] = {"0p1/2_1/2-_1", "0p1/2_1/2-_2", "0p3/2_3/2-_1", "0p3/2_3/2-_2",
  	"0p3/2_3/2-_3", "0d5/2_5/2+_1", "0d5/2_5/2+_2", "1s1/2_1/2+", "0d3/2_3/2+"};
  const double Ex[n] = {0., 3.663, 0.3739, 3.2042, 5.515, 2.526, 4.209, 0.8496, 5.195};
  const double c2s[n] = {1.6950, 0.0339, 0.2934, 0.3697, 2.5455, 0.0498, 0.0450, 0.0056, 0.0034};
#endif

#ifdef O19
  const char *basedir = "config/O/O19";
  const double A = 19., Sp = 17.0693, Sn = 3.9556;
  const double csExp = 14.3, dcsExp = 2.3;
  const int n = 9;
  const char *orbit[n] = {"0p1/2_2-_1", "0p1/2_2-_2", "0p1/2_3-", "0p3/2_1-",
  	"0p3/2_2-_1", "0p3/2_2-_2", "0p3/2_3-", "0d5/2_1+", "0d3/2_1+"};
  const double Ex[n] = {0.1147, 0.5874, 0.742, 0., 0.1147, 0.5874, 0.742, 2.6144, 2.6144};
  const double c2s[n] = {0.6849 ,0.0004, 0.9809, 0.0206, 0., 0.0387, 0.002, 0.0001, 0.0011};
#endif

#ifdef O20
  const char *basedir = "config/O/O20";
  const double A = 20., Sp = 19.3491, Sn = 7.6080;
  const double csExp = 21.68, dcsExp = 1.12;
  const int n = 4;
  const char *orbit[n] = {"0p1/2_1/2-", "0p3/2_3/2-_1", "0p3/2_3/2-_2", "1s1/2_1/2+"};
  const double Ex[n] = {0., 1.143, 2.132, 2.511};
  const double c2s[n] = {1.6661, 0.4974, 0.2763, 0.0073};
#endif

#ifdef O22
  const char *basedir = "config/O/O22";
  const double A = 22., Sp = 23.2379, Sn = 6.8503;
  const double csExp = 16.12, dcsExp = 1.1;
  const int n = 3;
  const char *orbit[n] = {"0p1/2_1/2-", "0p3/2_3/2-", "0d5/2_5/2+"};
  const double Ex[n] = {0., 1.160, 2.380};
  const double c2s[n] = {1.6315, 0.8734, 0.0213};
#endif


  PARA->SetBaseDir(basedir);
  // PARA->Clear(); // clean the previous text output of the base directory
  TARs rs(A, Sp, Sn, csExp, dcsExp); // the Rs calculator
  for(int i = 0; i < n; i++){ // loop over all the valence configurations
    char cfg[128]{}; rs.YamlName(orbit[i], cfg); // build the file name of the yaml
    PARA->Read(cfg); // read config file
    // PARA->SetValue("FNNParOpt", FNNPAR); // 0: Horiuchi; 1: LenziRay
    // PARA->SetValue("Ek", EK); // MeV/nucleon
    // PARA->SetValue("isPauli", PAULI); // 0: off; 1: on
    TAMOMDIS mom; // do the calculation
    rs.AddConfig(orbit[i], Ex[i], c2s[i], mom.Total()); // register the configuration
  } // end for over i
  rs.Rs(); // calcualte Rs

  return 0;
} // end the main function
